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ABSTRACT 

The masses and radii of extrasolar planets are key observables for understanding their 
interior, formation and evolution. While transit photometry and Doppler spectroscopy 
are used to measure the radii and masses respectively of planets relative to those of 
their host star, estimates for the true values of these quantities rely on theoretical 
models of the host star which are known to suffer from systematic differences with ob¬ 
servations. When a system is composed of more than two bodies, extra information is 
contained in the transit photometry and radial velocity data. Velocity information (h- 
nite speed-of-light, Doppler) is needed to break the Newtonian MR~^ degeneracy. We 
performed a photodynamical modelling of the two-planet transiting system Kepler-117 
using all photometric and spectroscopic data available. We demonstrate how absolute 
masses and radii of single-star planetary systems can be obtained without resorting 
to stellar models. Limited by the precision of available radial velocities (38 ms“^), 
we achieve accuracies of 20 per cent in the radii and 70 per cent in the masses, while 
simulated 1 ms“^ precision radial velocities lower these to 1 per cent for the radii and 
2 per cent for the masses. Since transiting multi-planet systems are common, this tech¬ 
nique can be used to measure precisely the mass and radius of a large sample of stars 
and planets. We anticipate these measurements will become common when the TESS 
and PLATO mission provide high-precision light curves of a large sample of bright 
stars. These determinations will improve our knowledge about stars and planets, and 
provide strong constraints on theoretical models. 

Key words: planets and satellites: dynamical evolution and stability - planets and 
satellites: fundamental parameters - stars: fundamental parameters - Planetary sys¬ 
tems. 


1 INTRODUCTION 

The mass and radius of extrasolar planets are usually ob¬ 
tained relative to those of the host star. Generally, stellar 
evolution tracks are used for field stars to infer their phys- 
* E-mail: jose-manuel.almenara-villaaujf-grenoble.fr ical parameters (mass, radius, age) from the spectroscopic 

(JMA); Rodrigo.Diazaunige.ch (RED) parameters (Teff, logf/, and metallicity) obtained from the 
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modelling of a stellar spectrum. If the star is transited by a 


planet, it is possible to estimate the stellar density (Sozzetti 


et al. 20071, that can replace logg, typically the most un¬ 


certain of the three atmospheric parameters, as input for 
the stellar models. The typical error is 5 per cent in the 
mass and radius (Wright et al. but usually these er¬ 

rors do not take into account the systematic errors of the 
models that can be up to 10 per cent ( Boyajian et al.|2012 |. 
Asteroseismology provides stellar radii and masses with a 


typical precision of 3 and 7 per cent, respectively (Huber 
et al.||2013bl. But this determination uses scaling relations 


(Ulrich 1986 Kjeldsen fc Bedding||1995 1, that require cali¬ 
bration using stellar models and also depends on the stellar 
effective temperature. 

The Gaia satellite ( Perryman et al.]|2001 1 will provide 
the stellar radius of a large number of stars with a precision 
of around 3 per cent, but the measurement will depend on 
the stellar effective temperature, the extinction, the bolo- 
metric correction, and the distance to the star. The best 
source of direct empirical determination of stellar masses 


and radii are double-lined eclipsing binaries (Torres, Ander- 
|sen fc Gimenez||2010[ ). Equivalent measurements in plane¬ 
tary systems imply the detection of the planet radial veloc¬ 
ity ( Snellen et al.|[2bl0 |, not achievable for a large number 
of planets at present. 

While the orbital elements of single-planet systems are 
constant, they are functions of time when more than one 
planet is present. In particular, planet-planet interactions 
perturb the timing of transits so that they are no longer 
strictly periodic. When detectable, the associated transit 
timing variations or TTVs and transit shape variations con¬ 
tain valuable information about the system parameters (e.g. 


Ragozzine & Holman 20101, especially planet-to-star mass 


ratios which are not available for single-planet systems. 

When photometric timing data alone is available, it 
is not possible to deduce absolute masses and radii using 
purely Newtonian modelling because of the latter’s inherent 
MR~^ degenerac3Q This degeneracy can be broken, how¬ 


ever, when velocity information is available (Agol et al. 2005 
Montet fc Johnson|2013 1, either in the form of variations in 


the time-of-arrival of light signals due to barycentric motion 
in a triple or higher-order system, or in the form of radial 
velocity measurements (or both). This information provides 
an absolute scale of the system. While the timing effect due 
to the finite speed of light has been used to directly mea¬ 
sure the masses and radii of the three stars in the triple 
system KOI-126 ( [Carter et al.|[^ll[ ), it is currently unde¬ 
tectable in single and binary stars hosting planets, being at 
the level of 1 s or less. However, supplementing photometric 
data with radial velocity measurements has allowed for the 
direct measurement of the masses and radii of both the stars 
and planets in the circumbinary planetary systems Kepler- 
16 ( [Doyle et al.|20li| , Kepler-34, and Kepler-35 ( [Welsh et al.| 


^ In other words, the model is invariant to scaling the lengths by 
a factor and the masses by the same factor at cubic exponent. 


20121. Note that the fractional uncertainties in the masses 


and radii of the three stars in KOT126 reduced from 10 and 
3 per cent, respectively, using photometry alone, to 3 and 
0.5 per cent when low-resolution radial velocity information 
was included. 

While the success of the studies described above is in 
large part due to the exquisite photometric precision of the 
Kepler telescope ( Borucki et al.|2010 |, an extremely impor¬ 
tant aspect was the use of photometricdynamical (or photo- 
dynamical) modelling. 

To date, most studies have proceeded in a different man¬ 
ner: the individual transits are fitted separately to obtain 
the times of mid-transit, often fixing the remaining transit 
parameters and thus effectively averaging over (and there¬ 
fore discarding) valuable information contained in the whole 
transit. The obtained transit times are used to obtain a mean 
ephemeris, the departures from which are finally fitted using 
an Ai-body integrator. This allows computing the N-body 
model at a much sparser resolution that the one needed to 
model each single photometric observation, reducing drasti¬ 
cally the required computing time. In contrast, the photody- 
namical approach models the whole light curve consistently 
by coupling N-body integrations with a model of the flux 
variations due to transits and occultations. 

The aim of this paper is to demonstrate that the photo- 
dynamical approach can be successfully employed to deter¬ 
mine accurate masses and radii of both the star and planets 
in single-stax planetary systems independently of any stellar 
models when high-precision radial velocity da ta are avail- 
ablepj Applying this ap proach to Kepler-117 ( Rowe et al. 
2014 [Bruno et al.|2015l, a system composed of an F9-type 


star, a 0.7 Rj planet in an 18.8 d orbit (planet b), and a 


1.1 Rj companion in a 50.8 d orbit (planet c). Bruno et al. 


(20151 detected TTVs in this system, which allowed them to 
measure the masses of the system planets using additional 
information on the stellar parameters based on spectroscopic 
measurements. 

As the light time travel effect is negligible in this system, 
we use the SOPHIE radial velocities reported by [Bruno et al.| 


(20151 to break the degeneracy mentioned above. In Sect. 2 
we describe the data employed in the analysis. In Sect. 3 
we present the details of the photodynamical modelling. In 
Sect. [^ we present our results, and finally in Sect. we 
discuss them and their importance in the framework of up¬ 
coming space missions like TESS and PLATO. 


2 DATA 

Kepler observed 67 transits of the inner planet and 29 tran¬ 
sits of the outer one between 2009 May and 2013 May. The 

^ Note that the photodynamical model has been used before 
to study multi-planetary systems with single host stars but the 
stellar parameters have been constrained using asteroseismology 
( [Carter et al.|[2012[ [Huber et al.[[T013a[ l, thus again relying on 
stellar models. 
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Kepler light curves of all Quarters (Q1 - Q17) were retrieved 
from the Mikulski Archive for Space Telescopes (MAST) 
archiv^We preferred short-cadence data (about one point 
per minute; quarters Q4-Q17) whenever available. Only the 
data around transits was modelled, after normalizing it with 
a parabola, and to account for the integration time of long 
cadence data (about one point every 30 min; quarters Ql- 
Q3), the model light curve was oversampled by a factor of 
20 and then binned back to match the cadence of the data 
points. We use the ‘SAP’ light curve, which we corrected for 
the flux contamination using the ‘CROWDSAP’ value in the 
hts hie header. We do not detect transit depth differences 
among different seasons. Therefore, the given contamination 
values seem self-consistent. 

The SOPHIE (Perruchot et al. 2008 Bouchy et al. 20131 
radial velocity observations of Kepler-117 are described in 


Bruno et al. 120151. Briehy, 14 radial velocity measurements 


at a resolution power of 40 000 were obtained between 2012 
July and 2013 November. 


3 PHOTODYNAMICAL MODEL 


The photodynamical model describes the light curve and 
radial velocity data at any moment in time accounting for 
the dynamic interactions via an Y-body simulation. The 
model parameters are the stellar mass and radius, the co- 


efhcients of a quadratic limb-darkening law (Manduca, Bell 
&: Gustafsson|1977[ ), and the planetary mass, planet-to-star 
radius ratio, and the orbital parameters (a, e, i, uj, n, and 
M; see Table at a hxed reference time (tref) for each 
orbiting planetary companion. The system was integrated 
over the span of Kepler and SOPHIE observations with a 
time step of 0.04 d, which produces a maximum error of 
2 ppm in the interpolated light curve model. In this way, 
the instantaneous system parameters, including the sky- 
projected planet-star separation, are known at each step. 
This, together with the planet-to-star-radius ratio and the 
limb darkening coefficients, defines the instantaneous light 
curve model ( [Mandel fc Agol]|2002[ ). Additionally, the line- 
of-sight projected star velocity is computed and compared 
to the SOPHIE radial velocities. As Y-body integrator, we 
employed the Bulirsch-Stoer algorithm implemented in the 
MERCURY code ( Chambers|1999 ). 

To sample from the posterior distributions of the pa¬ 
rameter models we used the Monte Carlo Markov Chain 


(MCMC) code implemented in the FASTIS package ( Dfaz 
et al.|20l4 1. We use the Huber et al. (2013aI parametrization 
to minimize the correlation between the model parameters, 
which can reduce the efficiency of the MCMC algorithm, 
and impede a proper exploration of the entire joint poste¬ 
rior distribution (see Table [^. 

We included a number of additional parameters in the 
model: a radial velocity linear drift, a global light curve 


^ http://archive.stsci.edu/index.html. 
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Figure 1. Two possible configurations for orbital inclinations. 
Stellar crossings paths seeing by Kepler of 1000 random MCMC 
steps from the photodynamical model fit are plotted. Planet c: 
red (inclination planet c > 90), and blue (inclination planet c < 
90), planet b: grey. Stellar disc in light grey. 


normalization factor for the short- and long-cadence Kepler 
data, and a multiplicative jitter parameter for each data set. 
The longitude of the ascending node at tref of only one of 
the two planets is explored, while the other is kept fixed 
(this is equivalent to fitting the difference of longitudes of 
the ascending nodes at tref). For a spherical star, the model 
does not depend on the values of the individual longitudes 
of the ascending node. 

Finally, we considered the two possible configurations 
for orbital inclinations: both planets transit the same or dif¬ 
ferent stellar hemispheres (Fig. [^. We set the planet b to 
transit one hemisphere {ib < 90°) and leave planet c free 
to be on any hemisphere (0° < ic < 180°). Uniform pri¬ 
ors were used for all the parameters. The starting point of 
the MCMC algorithm in parameter space was the previous 
solution found for this system ( Bruno et al.|[20T5 l. We ran 
40 chains of 50 000 steps each. The chains were thinned us¬ 
ing their autocorrelation length and merged together for a 
total of 3817 independent samples from the posterior. The 
mode and the posterior 68.3 per cent credible interval of the 
system parameters are given in Table 


4 RESULTS 

Figs |A1| and |A2| show the photodynamical transit model. 
Note that our model naturally reproduces simultaneous 
transits, as in transit no. 15 of planet c (no. 40 of planet b). 
In Fig. we show the SOPHIE radial velocities and the 
model 1, 2, and 3 a credible intervals for each time. Fig. |A3| 
shows the parameter posterior distributions and correla¬ 
tions. 
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Table 1. Model parameters. Posterior mode and 68.3 % credible intervals. The orbital elements have the origin at the star (Asteroidal 
parameters in the MERCURY code) and are given for the reference time = 2454967.63 BJDxdb- 


Parameter 


Mode and 68.3 % credible interval 

Stellar mass, M* [Mq] 


0 

^•^^-0.20 

[0; 1.557]“ 

Stellar radius, 


1.12 ± 0.23 

[0.490; 1.884]“ 

Stellar density, p* [p©]^ 


0.2886 ± 0.0065 


Surface gravity, log g [cgs] 


3.9737;5f 


Linear limb darkening coefficient, 


0.420 ± 0.026 


Quadratic limb darkening coefficient, 


0.130 ± 0.042 


Systemic velocity (at BJDxdb 2,456,355) 

, 70 [kms **] 

-12.9506 ± 0.0097 


Linear radial velocity drift, 71 [m s~^ yr“ 

-iji, 

3 ± 20 




Kepler-117b 

Kepler-117c 

Semi-major axis, a [au] 


0.101 ± 0.023 

0.197 ± 0.044 

Eccentricity, e 


^•^^^^'-0.0018 

0.03085lO;°°«4 

Inclination, i [°]^ 


88.667 ± 0.042 

89.644 ± 0.043“ 

Argument of pericentre, u [°] 


256.6^1-7 

300.i1®'5 

Longitude of the ascending node, n [°]^ 


181.18 ± 0.20 

180 (fixed) 

Mean anomaly, M [°] 


340.51';® 

141.0114;®® 

Radius ratio, Rp/R^,^ 


0.04689 ± 0.00014 

0.07105 ± 0.00017 

Planet mass, Mp [Mj] 


n 031 +0-031 

U.Ucii_o 014 

u.DO_q 32 

Planet radius, Rp[Rj] 


0.51 ± 0.12 

0.77 ± 0.18 

Planet density, pp [g cm“^] 


n Q01 q-|-0.0055 

U.o4i5_o.oi5 

1 '77'q”I"0-035 

1- ^ '^-0.037 

Planet surface gravity, logpp [cgs] 


9 c:Qec+0-059 

Z.ooo_ox2 

3.467lo;?2’’ 

ai [BJDtdb-2,450,000]''’*' 


4978.81123 ± 0.00058 

4968.63162 ± 0.00038 

02 


18.774480+© ^ '“”4 

— 1.4 X IQ--^ 

50.778307] X I0-4 

Kepler long-cadence jitter^ 


1.069 ± 0.036 


Kepler short-cadence jitter^ 


0.9962 ± 0.0024 


SOPR/F jitter*’ 


0-98la?4 


, b "1" , c ^ 

Af* 


0.00165031®'] ^ '^“7 


A-'-fp,b b 

A^PjC 


0 05151+^'® ^ '0“'* 
U.U5iai_^ 3 X 10-3 


ebcostub - ^ecCOSLOc *’ 


U.U^ZO_0 0033 


Cb cos a;b + — Cc cos cjc ^ 


0 Qi 70+0-0063 

U.Ui/o_0 0025 


Cb sin ujjj — sin cjc “ 


0.0006 ± 0.0019 


Cb sin a;b + — Cc sin cjc ^ 

“b 


-0.1034i1]2 X 10-3 


^ 99 per cent Highest Density Interval (HDI). 

^ MCMC jump parameter. 

reflected with respect to 2 = 90°, the supplementary i 
^ ai = tj-ef ~ ^ — E e sin E) with E = 2 arctan 

the orbital periods of the planets, that moreover is not 
orbital parameters at tj-ef should be used in a A^-body in 
Mq = 1.98842 Xl03° kg, Rq= 6.95508 xlO® m, Mj = 

mgle is equally probable. 


\ \ l 7^ tan \ 

y y 1 + e L 

constant, oli 
tegration. 

1.89852 xlO 

2(2 “2 - ; "2 IS not representative of 

and 012 should not be used to predict transit times. Instead, the 

kg, Rj = 7.1492 xlO'^ m 


The model parameters are in agreement with the val¬ 


ues reported by Bruno et al. (20151, but the uncertainties 


are significantly reduced: the uncertainties of the unitless 
parameters are between 2 and 28 times smaller. As we dis¬ 
cuss below, this is one of the advantages of employing the 
full photodynamical model instead of computing the central 
times of transit separately. 

With our model the information in the TTVs and the 
transit shape is fully exploited, which explains the improve¬ 
ment with respect to the|Bruno et al.j(2015l parameter de- 
“ " l|4p ^ “ 


termination. Figs[^and|^show the evolution of transit times 


and shapes, respectively. As can be seen, changes in tran¬ 
sit shape are clearly detected for planet b. On the other 
hand, the transit shape variations predicted by the model 
for planet c are too small to be detectable with the available 
data. The transit duration of the inner planet increases by 
14.7±1.7 min during the timespan of the if epZer observations 
(Fig.[5|. The duration and depth of the transits of planet b 
are increasing as the transits become progressively central 
(Fig. |A7| ). Eventually the transit paths along the stellar disc 
will cross, and mutual eclipses between planets could occur. 
This effect can aid to distinguish between the two possible 
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Figure 2. Radial velocities of Kepler-117 observed by SOPHIE a,nd photodynamical model fit. The 68.3, 95.5, and 99.7 per cent Bayesian 
confidence intervals are plotted in different intensities of grey. 



Time [minutes] 


Figure 3. TTV of planet b. Transit ingress and egress observed 
by Kepler, with relative flux colour coded (the so-called river or 
lava plot). The times are measured from the expected transit time 
for a linear ephemeris. The transits are sorted according to the 
orbital phase of the 1:1 resonance period. 

configurations of the orbital inclinations, either both plan¬ 
ets in the same or different halves of the stellar disc. With 
the available data both configurations are equally likely. An¬ 
other possibility would be to precisely measure the transit 
duration of planet b in a few years from now. The model 
predicts that by 2020 the transit duration of planet b will 
differ by ~ 3.6 min between both configurations. 

The central times of the transit obtained from the pho¬ 
todynamical model exhibit a 29 min amplitude variation 
with respect to a linear ephemeri^for the inner planet. We 
also detect clear TTVs of the outer planet of 3 min ampli¬ 
tude (Fig. [^. 

The linear ephemeris named trough the paper refers to a linear 
fit to the mid-transit times (of the transits observed by Kepler) 
obtained with the photodynamical model fit (Sect. [^. The lin¬ 
ear ephemeris for planet b is: BJDtdb = 2 454 978.8214(12) 
-h 18.795931(27) x Epoch, and for planet c: BJDtdb = 2 454 
968.63220(25) + 50.790374(15) x Epoch, the errors of the last 
digit are indicated in parenthesis. 


Epoch Planet b 

0 5 10 15 20 25 30 35 40 45 50 55 60 65 70 75 



Time [days] 

Figure 4. Transit shape evolution of planet b during Ke¬ 
pler observations. Upper panel: Phase-folded Kepler photome¬ 
try and maximum-likelihood transit model corrected for TTV for 
planet b. The transit models of all epochs are plotted but are 
hardly distinguishable. Lower panel: difference between the tran¬ 
sits of different epochs (colour coded) and the model of the first 
transit. The effect is observed in the data points by combining 
one third of the transits (black lines and circles with 3 different 
colours). Note that the difference between different epochs reaches 
600 ppm. 
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Figure 5. TDVs of planet b and c during Kepler observations. The two possible solutions for the inclination of planet c are plotted in 
different colours, ic >90° is shown in red and ic <90° is shown in blue (slightly shifted in epoch for clarity). The duration change during 
Kepler observations (last - first transit) is: planet b = 14.6 di 1.7 {ic > 90), 14.8 di 1.7 {ic < 90) minutes; planet c = 25 di 76 {ic > 90), 
—53 di 77 {ic < 90) seconds. The mean error in the transit times of planet c during Kepler observations is 20 s, so it is not possible to 
differentiate between the two configurations. 
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Figure 6. TTVs of planet b (upper panels) and planet c (lower panels) obtained with the photodynamical model fit. From left to right: 
TTVs (calculated respect to a linear ephemeris) against mid-transit time. Lomb-Scargle periodogram ( [Press &: Rybic^|1989| l of the 
TTVs. The peaks corresponding to the modulation periods of different resonances are noted (Pi;i = 29.83di0.23 d, P2:i = 72.3dil.3 d, 
= 170.9di7.4 d, and P4;i = 39.15di0.39 d). The remaining peaks correspond to aliases of the sampling frequencies (see [Roberts, j 
[Lehar Sz Dreher|1987[ |. TTVs folded at period of the 1 :1 resonance (29.83 d). 


With a period ratio of around 2.7, the Kepler-117 sys¬ 
tem is many resonance widths away from both the 2:1 and 
3:1 resonances, the widths of whic h are around 0. 06 and 
0.02, respectively, for this systenf] ([Mardliiig 20131. None 


the less, these two strong harmonics of the disturbing func¬ 
tion are responsible for significant perturbations to the or¬ 
bital elements (Fig. |A6| ), which in turn, manifest themselves 
in detectable TTVs, as well as significant Fourier power in 
the associated circulation periods of the harmonic angles. 


® The width of a resonance, Act, is defined to be such that the 
associated harmonic angle librates (rather than circulates) when 
CT — Act < ct < ct -|- Act, where ct is the period ratio. Physically, a 
librating harmonic angle allows for substantial energy and angular 


momentum to be exchanged between the orbits over many orbital 
periods, as is evidenced by the significant TTVs one associates 
with resonant systems. 
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For a system with orbital periods Ph and Pc and period ra¬ 
tio a — Pc/Pb, the circulation period of an n': n harmonic 
is Pn':n = Pc/\na — n'|, which for the 1:1, 2:1, 3:1, and 4:1 
harmonics of the Kepler-117 system are 29.83 d, 72.27 d, 
170.86 d, and 39.15 d respectivelj[^ matching the periods de¬ 
tected in the Lomb-Scargle analysis of the TTVs (Fig.[^. In 
spite of being far from resonance, the TTVs of systems like 
Kepler-117 are rich with information about their architec¬ 
ture and mass distribution; this will be explored elsewhere 
(Mardling, in preparation). 

We emphasize that the observed TTVs and transit du¬ 
ration variations (TDVs) are based on transit times and du¬ 
rations that are obtained as a by-product of the photody- 
namical code (obtained, respectively, as the mean and the 
difference of the first and fourth contacts computed from the 
sky-projected planet-star separation). These measurements 
are based on the data but are assisted by the assumption 
that they are produced by the gravitational interactions be¬ 
tween the system bodies. Thus the achieved precision in the 
determination of transit times are improved by a factor of 
4 with respect to a classical determination using only the 
individual transit light curves. 

Finally, the absolute masses and radii of the star and 
both orbiting planets are constrained without using stel¬ 
lar evolutionary models, mainly because of the effect in the 
transit times. The precision achieved for the star radius and 
mass is 20, and 70 per cent, respectively, limited by the 


available radial velocity precisiorj^ 


Bruno et al. 


(20151 re¬ 


sort to stellar models and obtain a much better precision. 


The discrepancy between Bruno et al. 120151 and the values 


reported here is signihcant to 95 per cent. As the available 
radial velocity data are not selective enough the posterior 
distributions are influenced by the chosen priors. In this case 

- as the stellar density and radius prior distributions are flat 

- the mass prior probability increases towards lower values, 
which explains, at least partially, the bias for less massive 
stars. The radius is in turn affected via the correlation with 
the stellar density. On the other hand, the parameters that 
are precisely constrained by the available data are in agree¬ 


ment with Bruno et al. 120151. In the next section we discuss 
the potential of this technique for stars for which 1 ms“^ 
precision is achievable. In this case, we will show that the 
mass posterior distribution is completely dominated by the 
synthetic data, and that the obtained mass is in agreement 
with the simulated value. 


° using as periods the median time between consecutive transits, 
see Fig[X4| 

^ The mean radial velocity error is 38 ms“^, to be compared 
with the radial velocity amplitudes reported by |Bruno et al.|2015| 
~7 ms“^ for planet b and ~90 ms“^ due to planet c. 


4.1 Adding high-precision synthetic radial 
velocities 

To probe the capability of the method to characterize targets 
coming from next generation space-based transit surveys, 
we replaced the SOPHIE observations by simulated radial 
velocity data with a precision of 1 ms“^(see Fig. [^. The 
parameters used for the simulation are listed in Table 
The stellar mass was chosen close to the one presented by 


Bruno et al. (20151, but this is irrelevant for the comparison 


below. 

The photodynamical modelling was repeated identically 
and the results are shown in Table With this precision in 
the radial velocities, the stellar radius and mass are mea¬ 
sured with a precision of 1 and 2 per cent, respectively. As 
the planet-to-star mass and radius ratios are known to a 
precision better than 2 per cent, independently of the radial 
velocities, the absolute radii and masses of the planets are 
also known to 1 and 2 per cent precision. 


5 DISCUSSION 

We have presented the analysis of Kepler-117 modelling the 
dynamical evolution of the system during the timespan of 
the Kepler observations using an A-body simulation. Usu¬ 
ally, the works in the literature studying the dynamical in¬ 
teractions in multiplanet systems compute first the transit 
times assuming a fixed transit shape at each epoch and then 
model the deviation of the thus measured transit times from 
a linear ephemeris (the TTVs). However, a number of ad¬ 
vantages exist in employing a complete dynamical model of 
the system over this two-step method. 

(a) Transit times are obtained consistently with the in¬ 
teractions and dynamical evolution of the system. This leads 
to a much better precision in the transit times. We refer to 
this determination as gravitationally-assisted transit times 
(see Fig. [^. In the case studied here, the transit times are 
determined with uncertainties four times smaller than the 


ones reported in Bruno et al. (20151. 


(b) As a consequence of (a), the precision in transit pa¬ 
rameters is improved. Including the shape changes in the 
model leads to an improved determination of the transit pa¬ 
rameters and the derived quantities. With respect to the 
analysis of [Bruno et al.| ( |2015[ ) the stellar density is deter¬ 
mined with twice the precision, the densities of planets b 
and c are known 11 and 4 times (respectively) better, and 
the difference of the longitude of the ascending node has an 
uncertainty 28 times smaller. 

(c) The masses and radii of the objects are obtained in¬ 
dependently of stellar models if the system scale can be 
determined. To do this, we have resorted to the SOPHIE 
radial velocities, but detecting the light time travel effect, 
which is negligible for Kepler-117 is another possibility. The 
independence from stellar evolutionary models make these 
determinations as valuable as those obtained in double-lined 
eclipsing binaries. 
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Table 2. Model parameters for the simulated run with lms“^-precision radial velocities. Simulated value, posterior median and 68.3 % 
credible intervals. The orbital elements have the origin at the star (Asteroidal parameters in the MERCURY code) and are given for the 
reference time = 2454967.63 BJDtdb- 


Parameter 

Simulated value 

Median and 68.3 % credible interval 

Stellar mass, M* [Mq] 

1.105 

1.119 ± 0.017 [1.0736, 1.1575]“ 

Stellar radius, [-R©]^ 


1.561 ± 0.012 [1.5321, 1.5865]“ 

Stellar density, p* [p©]^ 


0.2948 ± 0.0048 

Surface gravity, log g [cgs] 


4.1010 ± 0.0057 

Linear limb darkening coefficient, 


U.4U4_Q Q;J^4 

Quadratic limb darkening coefficient, 


U.105_o 044 

Systemic velocity (at BJDtdb 2,456,355), 70 [kms“^]^ 

-12.94774 

-12.94779 ± 0.00029 

Linear radial velocity drift, 71 [m s~^ yr~^]^ 

-12.62 

-12.25 ± 0.59 


Kepler-117b 


Semi-major axis, a [au] 

0.14302 

0.14362 ± 0.00074 

Eccentricity, e 

0.05263 

0.05251 ± 0.00085 

Inclination, i [°]^ 

88.700 

88.708 ± 0.040 

Argument of pericentre, u [°] 

256.5 

256.3 ± 1.9 

Longitude of the ascending node, n [°]^ 

181.31 

181.22 ± 0.18 

Mean anomaly, M [°1 

340.5 

340.7 ± 2.2 

Radius ratio, 


0.04679 ± 0.00012 

Planet mass, Mp [Mj] 

0.0932 

0.0946 ± 0.0016 

Planet radius, Rp[Rj] 


0.7103 ± 0.0063 

Planet density, pp [g cm“^] 


0.3275 ± 0.0073 

Planet surface gravity, logpp [cgs] 


2.6671 ± 0.0073 

ai [BJDtdb-2,450,000]^’^ 

4978.81123 

4978.81127 ± 0.00048 

CK2 

18.774380 

18.774386 ± 4.3 x lO"® 


Kepler-117c 


Semi-major axis, a [au] 

0.2776 

0.2788 ± 0.0014 

Eccentricity, e 

0.03099 

0.03121 ± 0.00073 

Inclination, i [°]^ 

89.685 

89.685 ± 0.034 

Argument of pericentre, u [°] 

300.2 

300.0 ± 1.4 

Longitude of the ascending node, n [°] 

180 

180 

Mean anomaly, M [°1 

140.9 

141.1 ± 1.5 

Radius ratio, R^i/R* 


0.07095 ± 0.00014 

Planet mass, Mp [Mj] 

1.828 

1.849 ± 0.027 

Planet radius, Rp[Rj] 


1.0771 ± 0.0091 

Planet density, pp [g cm“^] 


1.838 ± 0.035 

Planet surface gravity, logpp [cgs] 


3.5973 ± 0.0067 

ai [BJDtdb-2,450,000]'*’*’ 

4968.63151 

4968.63152 ± 0.00031 

«2 [d]'*’*’ 

50.778361 

50.778345 ± 8.0 x 10“® 

Kepler long-cadence jitter^ 


1.067 ± 0.039 

Kepler short-cadence jitter^ 


0.9963 ± 0.0024 

Radial velocity jitter^ 


1.07 ± 0.25 

M* 

0.0016591 

0.0016584 ± 2.9 x 10“® 

Mp b 5 

Mp c 

0.05099 

0.05112 ± 0.00030 

Cb cos U)-^ — cos U)c ^ 

-0.0425 

-0.0425 ± 0.0017 

Cb cos i^b + oL ^ 

0.0180 

0.0178 ± 0.0028 

Cb sin ojiy — ^Ec sin cOc ^ 

0.0008 

0.0013 ± 0.0016 

eb sin tab -f ^ ec sin Wc 

-0.1032 

-0.1036 ± 0.0015 


In practice, the TTVs and the changes in the transit 
shape cannot be detected for all multi-planet systems. How¬ 
ever, we note that the TTVs of Kepler-117 were not consid¬ 


ered significant by other authors (Mazeh et al. 2013 Kip¬ 


ping et al. 20151. Other similar systems might have been 


overlooked as well. Even if computationally expensive, the 
photodynamical model, as the one described here, permits 


to fully exploit the observations of multi-transiting planets 
such as those obtained by the Kepler mission. 

The only assumptions of the model are the Newtonian 
Law of gravitation and the geometry of an opaque disc 
occulting a bright one that follows a limb-darkening law. 
However, we have identified five neglected effects: (a) the 
uncertainty in the contamination of the Kepler photomet- 
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Figure 7. Idem Fig.j^for the photodynamical model fitted with simulated radial velocities of 1 ms“^ precision. 
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Figure 8. Mid-transit time error of planet b (dots) and c (plus 
markers) until 2020 from the photodynamical model fit. The ver¬ 
tical dashed line mark the end of Kepler observations. Note that 
the precision in transit timing is minimum at the middle of Ke¬ 
pler observations. So, precision in transit timing depends on a 
combination of photometric precision, photometry sampling, and 
dynamics (through gravity). 
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ric mask (it is provided without error), (b) stellar activity, 
(c) relativistic effects, (d) the light-time effect, and (e) the 
non-sphericity of the objects. We have repeated the analy¬ 
sis adding a flux contamination (a) as a free parameter. We 
found an additional contamination factor of -0.7±1.5 per 
cent with respect to the one reported in the Kepler archive. 
The planet-to-star radius ratios seem to be the only param¬ 
eters affected, with distributions 3-4 times wider compared 
to the analysis with fixed contamination. In any case, it 
seems that the contamination of the photometric aperture 
can be directly measured in this kind of analysis. About the 
activity (b), there is a 0.1 per cent amplitude variability 
curve plausibly caused by stellar spots ( |Bruno| 
, that can affect the planet-to-star radius ra¬ 
tios and the radial velocity measurements. In general, stel¬ 
lar activity can have an effect on the parameters obtained 


m the light 


et al. 2015 


by the photodynamical modelling, but probably less than 
if transit times are measured individually. Here, instead, 
gravitationally-assisted transit times are constrained by the 
entire light curve. General relativity (c), and the finite speed 
of light (d) are negligible in this case. Nor do we expect a 
significant deviation from sphericity (e) as the star is rotat¬ 
ing slowly ( Bruno et al.|2015 l, and the planets are in distant 
orbits, so tidal effects are negligible. An additional source of 
error are other bodies in the system not taken into account. 
Their influence depends mainly on their mass. 

In our analysis, the precision achieved for the mass and 
radius with the available SOPHIE radial velocities is very 
poor, specially the mass. Current high-precision spectro¬ 
graphs like HARPS-N csLTi achieve ~4 ms”'^ precision radial 
velocities on this star. Future space transit search missions 
will focus on brighter targets for which radial velocity pre¬ 
cision better than 1 ms”^ will be achieved. Thus, we have 
simulated radial velocity measurements with a precision of 
1 ms“^ at the times of SOPHIE observations. We have re¬ 
peated the analysis with these data obtaining a precision 
of 1 per cent for the radii and 2 per cent for the masses 
(both stellar and planetary). As the ratios planet-star ra¬ 
dius and masses are determined very accurately (<2 per 
cent), the precision in the masses and radii of the star and 
planets is similar. This simulation shows the potential of 
this technique. We have shown that it is possible to achieve 
a precision comparable with the also direct empirical de¬ 


terminations in double-lined spectroscopic binaries (Torres, 
Andersen fc Gimenez|[20l0 K Besides, in binary systems the 


stars can be affected by interactions between components, 
especially in low-mass stars. 

When a planetary system exhibits detectable gravita¬ 
tional interactions, the orbital parameters of the system, 
the mean densities of all bodies the planet-to-star radius 
and mass ratios, and planet mass ratios are determined by 


The planet density can be written as pp = p* 



When dynamical interactions are detected all terms are known 
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the photometry alone, without resorting to stellar models 
or to further measurements. If the effect of the finite speed 
of light is not detected, there are two possibilities to obtain 
absolute masses and radii: (a) impose priors on the stellar 
radius or mass (e.g. based on spectroscopic measurements, 
asteroseismology or stellar radius estimated with Gaia), or 
(b) break the Newtonian degeneracy using additional mea¬ 
surements (e.g. radial velocities). In this paper we took the 
second approach. The shape of the radial velocity curve, ex¬ 
cept its amplitude, can be inferred from the analysis of the 
photometry. Therefore, in principle, only two radial velocity 
observations (ideally at amplitude extrema) are needed to 
obtain a measurement of the system scale, and determine 
the masses and radii of the system bodies independently of 
stellar evolution models. However, more than a couple of ra¬ 
dial velocity measurements is desirable, as the presence of 
drifts (in the case of Kepler-117 there is no detectable drift) 
or spectrograph systematics can bias the measurement. Stel¬ 
lar jitter due to spots can affect and limit the precision in 
radial velocity measurements, but its effect can be reduced 
observing with an infrared spectropolarimeter like SPIRou 
(Delfosse et al. 20131. Given that the full characterization 
of the system with radial velocities is not needed, this can 
reduce dramatically the amount of follow-up observations 
needed for future space missions like TESS and PLATO. 

In principle the same technique presented here is ap¬ 
plicable to all transiting multiplanet systems, limited by 
the amplitude of the radial velocities, the amplitude of the 
TTVs, the photometric precision and time sampling. In 
practice the main limitation comes from the Doppler preci¬ 
sion obtained on the relative faint transits host stars. In the 
future, missions like TESS ( [Ricker et al.|20l^ and PLATO 
I Rauer et al.|2014 1 will provide a large sample of multiplanet 
systems around bright stars. For quiet bright st ars a preci- 
sion better than 1 ms“^ is already achievable (Pepe et al. 


2011). However, multiplanet systems are mostly composed 


of small planets (Mullally et al. 2015). Therefore, even for 
bright transiting hosts, radial velocities will limit the appli¬ 
cability of the method or the precision at which absolutes 
masses and radii are determined. Other limiting factor is the 
time sampling of the light curve. The 25 s observing cadence 
of PLATO will be an improvement with respect to previous 
space missions, and furthermore the fast cameras, with 2.5 s 
cadence, can improve the transit timing precision for the 
brightest targets. 

In the context of PLATO, the photodynamical mod¬ 
elling presented here could be complementary to asteroseis¬ 
mology in determining the physical parameters of target 
stars, as not all stars with multiplanet systems will have 
detected pulsations. This mass and radius determination is 
only limited by the precision in photometry and radial veloc¬ 
ity measurements, as opposed to the ones determined using 


precisely and therefore the planet density can be obtained without 
further measurements (see Table [ijl. 


stellar models, that depend on the understanding of the in¬ 
volved physical processes taking place in the stars. 
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Figure Al. Transits of Kepler-117b observed by Kepler. Each panel shows the photometric time series of the transits observed by Kepler 
and modelled in this work. In the lower part of each panel the residuals after subtracting the model to the observed data are shown. 
Each panel is labelled with the transit epoch. Each transit is centred relative to a linear ephemeris. 
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Figure A2. Idem Fig. |Al|but for planet c. 
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Figure A3. Two-parameter joint posterior distributions of the MCMC model parameters. The 39.3, 86.5, and 98.9 per cent two-variable 
joint confidence regions (in the case of a Gaussian posterior, these regions project on to the one-dimensional 1, 2, and 3 cr intervals) are 
denoted by three different grey levels. The histogram of each parameter is shown at the top of each column, except for the parameter on 
the last line that is shown at the end of the line. Units are the same as in Tabled 
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Figure A4. Time difference between consecutive transits, ‘period’, during A'ep/er observations for planet b (left) and planet c (right). 
The median value for planet b is 18.792 275 d, and for planet c is 50.790 557 d. 




Figure A5. Transit window of planet b (left), and planet b (right) for 500 years after Ti'ep/er observations. The 68.3, 95.5, and 99.7 per 
cent Bayesian confidence intervals are plotted in different intensities of red (inclination planet c > 90), and blue (inclination planet c < 
90). The solid line mark the median of the posterior distribution. T-axis positive means transit. Both planets will transit for the next 
500 years. 
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Figure A6. Evolution of the planets orbital parameters during Kepler observations from the photodynamical model fit (central body 
origin in the MERCURY code). The 68.3, 95.5, and 99.7 per cent Bayesian confidence intervals are plotted in different intensities of red 
(inclination planet c > 90), and blue (inclination planet c < 90). The solid line mark the median of the posterior distribution. Next to 
each plot, the Lomb-Scargle periodogram is shown. 
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The orbital parameters oscillate in the long term with a periodicity around 155-172 years. 

































































































































































